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ABSTRACT 

We study the effect of the non-Gaussianity induced by gravitational evolution upon the statistical 
properties of absorption in quasar (QSO) spectra. Using the generic hierarchical ansatz and the log- 
normal approximation we derive the analytical expressions for the one-point PDF as well as for 
the joint two-point probability distribution (2PDF) of transmitted fluxes in two neighbouring QSOs. 
These flux PDFs are constructed in 3D as well as in projection (i.e. in 2D). The PDFs are constructed 
by relating the lower-order moments, i.e. cumulants and cumulant correlators, of the fluxes to the 3D 
neutral hydrogen distribution which is, in turn, expressed as a function of the underlying dark matter 
distribution. The lower-order moments are next modelled using a generating function formalism in 
the context of a minimal tree-model for the higher-order correlation hierarchy. These different ap- 
proximations give nearly identical results for the range of redshifts probed, and we also find a very 
good agreement between our predictions and outputs of hydrodynamical simulations. The formal- 
ism developed here for the joint statistics of flux-decrements concerning two lines of sight can be 
extended to multiple lines of sight, which could be particularly important for the 3D reconstruction 
of the cosmic web from QSO spectra (e.g. in the BOSS survey). These statistics probe the under- 
lying projected neutral hydrogen field and are thus linked to "hot-spots" of absorption. The results 
for the PDF and the bias presented here use the same functional forms of scaling functions that have 
previously been employed for the modelling of other cosmological observation such as the Sunyaev- 
Zel'dovich effect. 



1 INTRODUCTION 

Ongoing Cosmic Microwave Background (CMB) experiments such as Planclfl the Atcama-Cosmology Telescop^ (ACT) and the South Pole Tele- 
scop^ (SPT) will pinpoint the cosmological parameters that describe the background geometry and dynamics of the Universe in an unprecedented 
detail. Along with very precise constraints on the structure of the Universe on the largest scales, smaller-scales observables will be crucial in order to 
further constrain the cosmological concordance model or find possible deviations from it. In particular, galaxy clustering (BOSS, 6dF, etc.) surveys 
and future weak lensing and clustering observations (e.g. EUCLID) could probe smaller scales and new redshif regimes. Spectroscopic surveys such 
as BOS:fl (and BIG BOSS) will also trace the large scale distribution of the baryonic matter in the Universe through the study of the flux distribution 
of the Lyman-a absorption systems in a very large number of quasar (QSO) spectra. 

The Lyman-a "forest", the many absorption features in QSO spectra, produced by intervening neutral hydrogen in the intergalactic medium 
(IGM) along the line-of-sight, is well known to be an important cosmological probe (for a recent review see Meiksin 2009). In the standard cos- 
mological paradigm, the IGM consists of mildly non-linear gas, making up the cosmic web, that traces the dark matter and is photo-heated by 



^ http://www.rssd.esa.int/Planck 

^ http://www.physics.princeton.edu/act/ 

^ http://pole.uchicago.edu/ 

http://cosmology.lbl.gov/BOSS/ 
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Figure 1. The PDF p{Fa) and bias b{Fa) of the flux Fa is plotted as a function of the flux Fa - The PDF (left-panel) and the bias (right-panel) both ai'e constructed 
using a lognormal model for the underlying mass distribution. The resulting PDFs for <5 are next transformed into the flux PDFs using the fluctuating Gunn-Peterson 
approximation Eq.@- We compare the results from the lognormal approximation (solid lines) ag ainst the one based on the Gaussian approximation (dashed lines). The 
values of A{z) and /3 are constructed using the functional fit in Eq.(9) given by IKim et alJl2Q07l) see text for more details. The bias changes signs at an intermediate 
flux value which depends on cosmology and hydrodynamical parameters that define the equation of state of the photo-ionized medium A{z) and /3. A fiducial value 
for the variance cr = 3 was assumed for this plot. 



a Ultra Violet (UV)-background. The Lyman-a forest is thus the main probe of the IGM and it has been shown to arise naturally in hierarchi- 
cal structure formation scenarios. Astrophys ical effects produced by feedback from ga laxies and/or AGNs do not seem to strongly affect the vast 
majority of the baryons in the cosmic web ^McDonald et alj|2005l : I Theuns et al ] l2002h . thereby this can be used as a dark matter tracer. The re- 
lation between the Lyman-Q forest flux and the underlying matter field is a nonlinear one and it is generally expected that statistics of Lyman-a 
are biased relative to the underlying dark matter distribution. The Lyman-a forest has been studied using variety of analytical techniques such 
as the Zel'dovich approximation, which is valid in the quasilinear regime and often used in modelling of the nonlinear gravitational clustering 
tooroshkevi ch & Shandaririi |l977l : IHuI. Gnedin & Zhanglll997l : lMcGilllll99(]l : iMatarrese & Mohavedl25o2l). In addit ion, the lognormal approxi- 
mation IColes & Jones 1991) is frequently used to model the statistics of Lyman-a forest (lBil ll993l : iGnedin & Hui| [r996: Bi & Davidson 199"^; 
Roy Choudhury. Padamanabhan & Srianand 200 ll : l\^el et al [2002^. Models based on the hierarchical or scaling ansatz (Balian & Schaeffeiil 19891 : 
Bernardeau & Schaeffei 1992. ,199 9) fo r higher-order correlation functions have also been investigated in order to model the statistics of Lyman-a 

forest jValageas. Schaeffer & Silkll 19991) . 

In addition to analytical modelling, hydrodynamical simulations have also played a very important role in this field (e.g. ICen et alj h994h : 
iGnedin & Huil ( Il998l) : ICroft et al.l jl998l . [l999): Meiksin & White (2001)) and support the simple analytical picture. Thus, analytical schemes, 
including the ones that we develop here, can be calibrated using numerical simulations and this in turn allows to explore a large parameter space 
efficiently. However, numerical simulations are required to resolve the Jeans scale of the photo-ionized warm IGM, and this requirement typically 
means small box sizes that sample larger scales modes rather poorly. Indeed several semi-numerical prescriptions, which are not entirely based 
on hydrodynamical simulations, have also been developed to model Lyman-a flux and recover the correlation function from observed data sets 
Isiosaretal. 2011). 

The two most commonly used approaches in Lyman-a studies are based either on decomposing the information encoded in the transmitted flux 
via Voigt profile fitting or treating the flux as a continuous field. 

In the first approach, the shapes and c lustering properties of absorption lines fitte d by Voigt profiles have been investigated in v ariety of stud- 
ies involving the temperature of the IGM jSchave et al.lll99^ ; iMcDonald et alj|200ll) . in order to con strain the reioniza t ion histor y ("Theuns et ah 
20021: Ih ui & Haiman 2003 ), to measure the matter power spectrum and cosmological parameters jCroft et al.l 1 19991: 1 Viel. Haenelt & Spring^ 
2004; McDonald, Selia k & Burlesj|2006l : [McDonald & Miralda-Escudelll999l : iRolinde et alj|2003l : ICoppolani et al.ll2006l : iGuimaraes et al.ll2007t 
Viel & Haehnellt.2006) . 



By using the second set of methods, s t atistical properties of the flux such as the mean flux level, flux PDF, flux power spectrum 
iViel, Haenelt & Springel|[2004l : IViel et a l' 2008: Seliak, Slosar, McDonald 2006) and flux bi-spectrum are typically employed to explore flux statis- 
tics. FOTecamgle, it has been shown that the mean-flux level can be used to constrain the ampli tude of intergalactic UV background (Tytler et al 
l2004lBolton etal. 2005) while the flux PDF jMcOuinn et al.ll2009l : Isolton. Oh. Furianettdl 20091) is sensitive to the thermal evolution of the IGM 
(see also 'Lidz et al. (2006)). The flux power sp ectrum on the ot her hand can be used to constrain the cosmological parameters and nature of dark 
matter ( Croft et al. 2002). The flux-bispectrum jViel et alj[2009l) contains useful information about the primordial, as well as gravity-induced (i.e. 
secondary) non-Gaussianity. The data typically used in these investigations consists mainly of two different sets of QSO spectra: the SDSS low 
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resolution low signal-to-noise spectra and UVESA'LT or HIRES/KECK high resolution spectra. The number of SDSS spectra is about a factor 
~ 200 larger than that of high resolution samples though the later probes the smaller scales with greater accuracy. 

More re cently, it has been argued that BOSS-like QSO spectroscopic surveys could detect Baryon Acoustic Oscillations ( BAO) signatures a t 
high redshift jMcDonald & EisensteinlllOOTh and a sample of QSO pairs can constrain the geometry of the high redshift universe jMcDonalj|2003h . 
Furthermore, analysis of coincident absorption lines in QSO pairs can also allow departures from the Hubble flow and non-gravitational effects to 
be measured ( Ranch et al., 2005,). 

The SDSS-III/BOSS survey aims at identifying and observing more than 160,000 QSO over ^ 10, 000 square degrees within a redshift range 
z = (2.15 — 3.5). This survey is primarily design to studying baryonic acoustic oscillations by performing a full 3D sampling of the matter density. 
Such studies will also provide an unprecedented opportunity to study the clustering statistics using projected Lyman-a flux decrements of QSOs. 

In this paper, we will consider the Lyman-a flux decrement in two dimensions (2D) which is related the projected density of neutral hydrogen. 
The statistical study of projected Lyman-a flux decrement can be performed using one- or two-point PDF or their lower order moments. The PDFs 
contai n information of cumulants or th eir correlators to an arbitrary o rder and can be constructed using well-established machinery of hierarchical 
ansatz jMunshi. Coles. Melotlll999ffllMunshi. Melott. Coles! 1999cl) . Several authors have recently st udied the lo wer-order c umulants of Lym an-a 
forests and cross-correlated them against weak lensing convergence as well as to the CMB sky ( Vallin oTto et alj |20Q9: Valen tino et al.ll201ll) . The 
results presented here are complementary to such studies as we take into account the lower order moments to an arbitrary order not just for the 
one-point cumulants but also for their two-point counterparts or cumulant correlators. The cumulant correl ators are the two-point analogues of 
one-point cumulants and are already in use in different areas of cosmology, e.g. in analysis of galaxy surveys jSzapudi & Szalavlll997h . The lowest 
in the two-point hierarchy is the two-point correlation function. In the context of Lyman-a studies, t he two -point correlation function has already 
been introduced in studies involving two neighbouring line of sights (Viel et al. 2002; D'Odorico et al ]|2002l) . Our study generalizes these results to 
probe non-Gaussian correlation functions involving multiple line of sight. We mo del the statistics of underly ing neutral hydrogen distribution using 
lognormal distribution as well as an extension of perturbation theory approach jValageas & Munshil 12004) in 3D. Predictions from these models 
are then tested using hydrodynamical simulations at three different redshifts {z = 2, 3, 4). Next we use these results to build and test statistical 
description of the projected flux distribution. The results presented here can be generalized to cross-correlation studies involving external data 
sets and Lyman-a flux distribution. The correlati on functions (CCs) are equivalent to their Fourier (harmonic) space counterpart, the multispectra, 
recently introduced by jMunshi & Heavensll2009l) . 

The plan of this paper is as follows. In ^we introduce the notations and define the relevant quantities such as the transmitted flux and its 
relation to the underlying density contrast. In ij3]w6 gi^s details of the simulations that were used in our study. In i|4]details of modelling of the 3D 
flux are presented. In fj5T|we derive the lower order statistics for the flux in terms of that of the underlying density contrast. In i ]5.2| we provide a 
very brief introduction to the hierarchical ansatz which we use to model the statistics of underlying density contrast. In i ]5.3l we provide derivation of 
the PDF and the bias associated with the flux distribution and finally [j6]is left for discussion of results. We also provide a brief appendix introducing 
the lognormal approximation as well as the hierarchical ansatz. 



2 NOTATION 

We will be using the following form of the Robertson- Walker line element for the background geometry of the universe: 

ds^ = -cdt^ + a(t)(dr'^ + d\{r){de^ + sin^ edfj)"^)) (1) 

The angular position on the surface of the sky is specified by the unit vector Cl — {9, 0). The scale factor of the Universe is given by a{t). We have 
denoted the comoving angular diameter distance by dA{r) where r denotes the comoving radial distance to redshift z; dA{r) — K~^^^ sin(K^''^r) 
for positive curvature, dA{r) = (— K)"^''^ sinh((— K)^'''^r) for negative curvature and r for the flat universe. For a present value of of Ho and Q,m 
we have K = (JIm — l)Ho. The Hubble constant is denoted by Hq. The comoving radial distance r{z) at a redshift z can be expressed through the 
following LOS integration, with z' playing the role of intermediate redshift along the line of sight: 

(2) 

-RosJ^mil + z'yi + !^k(1 + z'Y + S^A 

Throughout, in our calculation, we will adopt Ho = lOO/i kms^^Mpc^^ with h = 0.7 and as = 0.87. Here Q,m, f^A and f^K are dark-matter, 
vacuume-energy and curvature contribution to cosmic density JIm + ^a + = 1. For a flat Universe Ok = 0. 




3 SIMULATIONS 

We use simulations run with the parallel hydrodynamical (TreeSPH) code GADGET-2 based on the conservative 'entropy-formulation' of SPH 
iSpringel et al.|[200^ . These consist of a cosmological volume with periodic boundary conditions filled with an equal number of dark matter and 
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Figure 2. The PDF of the density contrast p{S) as a function of the density contrast 5. The left {z = 2), middle {z = 3) and right {z = 4) pan els con'espond the to 
differ ent redshifts z studied. Two different analytical models are displayed the lognormal approximation (long-dashed) and the model proposed by ^Valage as & Munshil 
(short- and long-dashed) denoted as VM.The solid lines correspond to the results from simulations. Notice that the analytical results and simulations for positive 
(over-dense) 5 values of the PDF agree extremely well for all three redshifts investigated in this work. The PDF was estimated using 128^ giid and can resolve PDF as 
low as 0(10^®). However the PDF becomes increasingly dominated by the presence or absence of rare large overdensities. 
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Figure 3. The PDF of density contrast p(5) as a function of the density contrast 5. The left (z = 2), middle (z = 3) and right (z = 4) panels correspond to different 
redshift z. Two different analytical models are displayed the lognormal approximation (long-dashed) and the model proposed by (Valageas & Munshi 2004) (short- 
and long-dashed) denoted as VM. The PDF is constructed on a grid of 128^. The lowest probability that can be estimated on this is grid is 0(10^*'). However higher 
5 tails of PDFs are increasingly dominated by the presence (or absence) of rare over (under)dense objects. 



gas particles. Radiative cooling and hea ting processes we r e fol lowed for a primordial mix of hydrogen and helium. We assumed a mean Ultraviolet 
Background similar to that proposed bv lHaardt & Madaj 11996) produced by quasars and galaxies as given by with helium heating rates multiplied 
by a factor 3.3 in order to better fit observational const raints on the tempe rature evolution of the IGM. This background gives a hydrogen ionisation 
rate ri2 ~ 1 at the redshifts 2 = 2 — 4 of interest here ("Bolto n et al.l200^ . The star formation criterion is a very simple, one that converts all the gas 
particles into stars whose temperature falls below K and whose density con trast is larger than 1000 (it has been shown that the star formation 
criterion has a negligible impact on flux statistics). More details can be found in jViel, Haenelt & SDringelll2004h . 

The cosmological reference model corresponds to a 'fiducial' ACDM Universe with parameters, at z = 0, ilm = 0.3, f^A = 0.7, ilb = 0.05, 
ria = 1, and Ho = 70 km s~^ Mpc~^ and erg — 0.85. The initial conditions are generated at z = 49 using as linear matter power the one 
extracted using the CAMB software. We have used 2 x 512'^ dark matter and gas particles in a volume of linear size 60 comoving Mpch~ 1. The 
gravitational softening was set to 4 kpc in comoving units for all particles. The mass per particle is 1.1 x 1O**M0 h"^ and 1.8 x IO^Mq h'^ 
for DM and gas particles, respectively. This simulation should have sufficient resolution to properly reproduce most of flux statistics at least at z < 3 
and marginally at z = 4: the flux power should be in fact converged with this setup and also (marginally) the flux PDF. However, since we are 
dealing with projected quantities, we believe that numerical convergence issues are less severe than for other small scale flux observables. 

We select 128^^ grid points in the x ~ y plane and we interpolate along lines of sights parallel to the 2— axis using a Cloud-In-Cell (CIC) 
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Figure 4. We show the projected transmitted Lyman— a flux (top panels) and projected gas density (bottom panels) in the simulated volume at z=2,3,4 (left, middle 
and right columns, respectively). The gas density is extracted at 128'' grid points using a CIC or cloud-in-cell interpolation scheme: the values shown in the bottom 
panels are the projected <5 along 128'^ lines of sight (along the z-axis). The growth of the cosmic structures can be appreciated from high to low redshift. In the upper 
panels we show the mean value {Fa) along the same 128^ lines of sight, the values at the three redshift are all normalized to the {Fa{z = 2)} value and the spectra 
are simulated by using the exact definition of the transmitted flux and not approximations. Since the mean flux is a strongly evolving function of redshift, the growth of 
cosmic structure in Lyman-a flux is more difficult to interpret than in the corresponding density slices and simulations and/or models are needed. In the upper panels, 
voids correspond to regions of high-transmissivity, while dense regions produce absorptions. At 2; = 2 the (Fa) values fluctuate by about 10% around the mean, while 
at z = 4 this value becomes 20% (since the mean transmitted flux is about two times smaller at this redshift). The simulation refers to a hydrodynamical run with 
2 X 512^ gas and DM particles in a 60'^ (Mpc h^^)'' periodic comoving volume. 



algorithm. The Lyman-a flux is also computed in redshift space along the same lines of sight by using the exact definition of transmitted flux and 
not the Fluctuating Gunn-Peterson Approximation. Peculiar velocities, neutral hydrogen fraction and gas temperature are calculated along the lines 
of sight and a Voigt profile is used to obtain the mock quasar spectra. The density and transmitted flux are then projected along the simulated lines 
of sight in the 2— direction. 



4 THE LYMAN-a FLUX AND THE IGM IN 3D 

The fluctuating Gunn-Peterson approximation jGunn & Petersoijl 19651) allows us to relate the transmitted flux Fa{z, Q.) along a line of sight at a 
direction Cl = {9, (j>) at a redshift z with the fluctuation in neutral hydrogen density contrast 5: 

Fa(z, Cl) = exp + S(z, Cl)fj ■ (3) 

A{z) = 0.0023(1 + zf-^'"; /3 = 2 - 0.7(7 - 1); 7 = 1-3. (4) 

The parameters A{z) and /3 are two redshift-dependent functions relating the fl ux fluctuat i ons to the dark matter overdensities. The parameter A is 
related to the level of mean transmitted flux which is accurately measured (e.g. lKim et al It also depends on inputs from baryonic physics 

such as the baryonic fraction, IGM temperature, photo-ionization rate as well as background cosmological parameters. The power-law index 7 of 
the IGM temp erature-density relation determimes the value of the parameter [] and is relatively independent of redshift z jHui & GnedinllT997l : 
iKim et alj[2007n (this is of course assuming fluctuations in temperature due to e.g. reionization play a sub-dominant role). The allowed range of 
values for 7 is 7 = (0.7 — 1.5) we will use 7 = 1.3 for our calculation, which is a good approximation to the actual value of the hydrodynamical 
simulation. 
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4.1 PDF and Bias in 3D 



We use the lognormal distribution to model the distribution of 5\ the resulting flux PDF is obtained via a change of variable: 



dS 



dF^ 



dS^ 



(1) 



dSo 



dF, 



(2) 



(5) 



In our notation 5(i) = 5{Q,i, Zi) and a similar notation is adopted for the fluxes along different lines of sight Fa'^ = Fa{^i, Zi). For 3D studies 
we have compared flux from neighbouring lines of sight at the same redshift. We have compared the PDFs using this approach for three different 
redshifts 2 = 2,3, 4. The PDFs for p{5) are constructed using lognormal prescriptions (Eg.l lBlb and Ea. iB4t ). 

The lognormal distribution is also used in the analysis based on the fitted column density distribution of absorption lines. In such studies the 
lognormal distribution is introduced to relate 5 and number density of neutral hydrogen nni- In relating the flux distribution and 5 on the other hand 
we assume that 5 traces the underlying dark matter distribution which can be described well by lognormal distribution. The lognormal distribution 
here is being used to model the effect of gravity induced non-linearity on PDFs of matter density contrast 5. A similar but somewhat different 
approach was used bv lViel et al.1 ( |2002|) who used lognormal approximation to map the density field to flux. In our case we directly map the PDF 
of the 5 to that of the flux using simply a change of variable Eq.(|4). The co nstruction of the PDF of 5 is done using lognormal approximation as 
well as using extensions of perturbative methods jValageas & Munshilliooj) . These methods give near identical results. The hierarchical ansatz (to 
be introduced later ; see Appendix-B) will be used to o for the construction of one-point and two-point PDFs. It is interesting to point out that it has 
also been used by jValageas. Schaeffer & SiuJi 1999 1) to predict the clustering Lyman-a absorbers. It was used to construct an unifying model from 
the clustering of low column-density clouds to the collapsed dense damped systems. The methods that we present here is complementary to these 
studies as we directly probe the statistics of the transmitted flux. 



4.2 Numerical Results in 3D 



In Figure ([2} we have plotted the PDF of p{5) as a function of 5. The PDF was computed by binning the data on a 128^ grid. This allows us to 
resolve the PDF down to O(10~^). Different panels show the PDF for three different redshifts. We also show the results from different analytical 
calculations. The two-point PDF or 2PDF denotes the joint PDF of 5 at two different points. Typically the 2PDF can be written as a product of two 
IPDFs at these two-point and additional correction terms. The corrective terms represents the contribution due to correlation. In general in the large 
separation limit when the two-point correlation is weak compared to the variance (i.e. < C12 = (<5(xi)5(x2)) is the two-point correlation 
function and ^2 is its volume average), the 2PDF takes the following form in 3D: 



p(5^,S2)d5id52 = p(5i)p(<52)[l + b(<5i)eJ2b(<52)]d<5id<52. 



(6) 



This equation describes the two-point PDF (or 2PDF) , 52 ) of the underlying density contrast 5 in terms of PDFs p{5) and bias b{5). An identical 
relation also holds for 5 above the Jeans scale. I n general the (diffe rential) bias is difficult to estimate. To increase signal-to-noise it is useful to 
compute the integrated bias beyond a threshold jBemardeaull 1 9941) : 



fe(> '5o) 



Is^ Is^ P(^^'^2)dSidS2 



1/2 



(7) 



For our numerical estimation of bias associated with the S we have used Eq.lO to compute the cumulative bias associated with 5. In Figure ([5) 
we have presented the numerical results for 6(> Sg) against the theoretical expectations from lognoraml distribution. As expected the lognormal 
distribution is quite accurate in prediction of the numerical results for over dense regions. We show results for three different redshift relevant to high 
redshift Lyman-a studies. Notice that the high 5 tail of the PDF is reproduced well by the lognormal approximation Eq.l |B5l >. The departure at lower 
5 is partly related to the smoothing introduced in hydrodynamical simulation during interpolation. Indeed we only use the leading order term ^12/^* 
in our calculation of bias from the lognormal distribution. These results are valid in the large separation limit. Additional corrective terms can be 
included if necessary. 

Using the transformation given in Eq.H, we can express the pdf p{Fa) and bias b{< Fa) in terms of the pdf p{So) and bias for the density 
bi> So): 



b{> So) 



Is"^ Is"^ P(^i^^2)dSidS2 



IL^ p{s)dsy 



- 1 



1/2 



[/,f>(F„)dF.] 



- 1 




(8) 



The results of such calculation for the PDF is shown in Figure There is no accurate model for bias in the highly nonlinear re gime. The lognormal 
appro ximation is accurate in the quasilinear regime and starts to fail in the highly nonlinear regime. The approach taken by jValageas & Munsha 
I2OO4I) can be extended to the construction of bias. The results of integrated bias for the flux as defined in Eq.([8]l and are plotted in Figure ([6). We 
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Figure 5. The PDF of the flux p(-Fa) as a function of the flux Fa. The left {z = 2), middle {z = 3) and right (z = 4) panels con'espond to different redshift z. 
Two different analytical models are displayed the lognormal approximation (long-dashed) and the Gaussian approximation (short and long-dashed). The flux PDF is 
constructed from the PDF of & using Eq.j4). The results from the simulation are also depicted (solid lines). 
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Figure 6. The bias of the flux 6(< Fa) as a function of the flux Fa- The left (z = 2), middle (z = 3) and right (z = 4) panels correspond to different redshift z. The 
solid lines coirespond to direct estimates from 3D flux maps. The dashed lines are from maps that are generated from the density maps. 



found that including a power law mapping (1 + 5) — >■ (1 + 5)^, where A is a constant of order unity, can fit the numerical results, if it is applied 
before the exponential transformation. 



5 PROJECTED FLUX DECREMENT 

We are interested in scales larger than the comoving Jeans length. At around redshift z — Z which is about l/i~^Mpc. The fractional flux 5 J- a {z, f2) 
in this limit can be related to the underlying dark matter density distribution 5, using Eq.l[4]l by the following expression: 

5Fa{z, Q) = {Fa{z, Cl) - {Fa{z, n)))/{Fa{z, Q)) = - A{z) l3{z)5{{l, z) ^ ^ A{z) {z)b{z)S {{l, z). (9) 

The relation assumes that IGM traces the dark matter distribution at larger scales (above Jeans length). We have incorporated an additional redshift 
dependent bias factor b{z) in relating Sigm{z) and 5{z) for generality. The results presented here can be generalised for arbitrary redshif t depen dent 
bias. In our calcu lations we will set b{z) — 1. The factors A{z) and /? have already been introduced above; see e.g. iHui & Gnedinl l ll997h and 
iMcDonaldl Jiooj) for a discussion on the temperature-density relation which is well described by a power law. 

The integrated flux decrement ST a towards a particular direction Cl can be expressed through a line-of-sight integration as a projected 
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density contrast: 

5J'^{Q.)= f " dr SF^{Cl,r) ^ f " drA{r)P{r)S{Cl,r)= f " dr uja{r) S{Cl,r). (10) 

The range of comoving distances probed by the Lyman-a spectrum extends from to the quasar distance r-q. We will define a parameter (5 J-"™'" = 
— '' Wa (r) dr that is the minimum value of SJ^a for a given ri and r,. The weight Wa {z) that acts as a kernel for projection can be expressed in 
terms of previously introduced functions w{z) — A{z)/3{z). Few comments are in order. The assumption of quasilinear theory ( e.g. the lognormal 
model) is expected to be valid above the Jeans length. This is where our results are expecte d to be valid (Eisenstein & Hu 1998). The relatio nship 
between the flux and the underlying mas s distribution has been verified in numerous studies jMcauiniTet a l. 2009; Bi & Davidson 1997; Cro ft et al.l 
l2002l ; IViel et alj2002l ; ISaitta et al.ll2008l) . It is also important to note that the simple linear relationship that is often used in deriving many analytical 
results in the context of Lyman-a studies can be modified in the presence of any non-gravitational process such as fluctuations in the level of 
ultraviolet radiation or temperature fluctuation. Such non-gravitational effects are harder to model analytically. Next we will use these expressions 
to model the lower order statistics of STc,{Cl). 



5.1 Lower order Cumulants and Cumulant Correlators in Projection 

We will use both cumulants and cumulants correlators in our studies. The one-point cumulants can be estimated from a single LOS where as 
neighbouring pairs of QSOs are required to probe the cumulant correlators. By construction, the normalised cumulant correlators are independent of 
the separation of quasar pairs. Increasing the number of line sights will help to increase the signal-to-noise. 

To compute the variance and other lower order moments or cumulants we start with the Fourier transform of the 3D density contrast 5 which 
we denote as 5(k): 

5Tc,{(l) ^ J^ " dr uj^{r) ^ exp(ir-fc|| + idAW^ia • kx)5(k, r). (11) 

Here 6^12 denotes the angle between the line of sight direction f2 and should be treated as a vector on the surface of the sky and the wave vector k, 
fc|| and k_L denote the components of k, parallel and perpendicular to the line of sight direction. In the small angle approximation one assu mes that 
Jkxl >> fc||. Using these definitions we can compute the projected variance (STa) in terms of the dark matter power spectrum Ps{k,r) jKaisej 
I1992I) : 

^^■^^(")^^^r'^^Sly/(£F^^(^'^) '''' 

Similarly the higher order moments of the field relate (SJ^^) to the tliree-dimensional multi-spectra of the underlying dark matter distribution Bp 
(Hui 1999; Munshi & Coles 1999a): 

Jr, 4 (27r)^ J (27r)2 \dA[r) dA(r) J 

The subscript c denotes the connected part of a diagram that represents the higher order correlation function; Limbers approximation ( lLimbeiiri954h 
is used to derived the above result. The subscript ~ 0) denotes a multiplicative Dirac's delta function fo[X]r=i We have used the 

compact notation B^*"' to denote the p-th order multispectra where B^'^' denotes the bispectrum and B^*' denotes the trispectrum, often denoted as 
Tg in the literature. The power spectrum correspond to B^'^' . We will use these results to show that it is possible to compute the complete probability 
distribution function (PDF) of SJ^a from the underlying dark matter PDF. The details of the analytical results presented here can be found in 
l lMunshi & JaiiJl200dl200lh . 

In addition to the cumulants their correlators are important as they are related to the two-point PDF (or 2PDF) and hence the bias associated 
with the higher flux regions 

i^mn.)snm>. ^ jyr-g^Srl-J^ . - . / i^^n'r" {j!^.-. ^)^^ , (», 

A special case of this equation corresponds to the expression for the two-point correlation function {5J-^{Cli)5J^a{^i)}- Tare generic. We will 
specialise them by assuming a particular form for the higher order multispectra to make further progress. Just as the normalised cumulants Sp 
we can define the normalised cumulant correlators Cpq = {5Ta{^i)SJ-i,{i^i))/{5J-a{i^i)SJ-a{^i))){5J-a{^))- We will use the normalised 
cumulants or the Sp parameters and the normalised cumulant correlators to construct the PDF and the bias associated with the flux decrements or 
ST a - We will denote the normalised cumulants and their correlators of density by Sp and Cpq. The corresponding quantities for the flux-decrements 
will be denotes using a superscript a i.e 5*^ and Cpg. 
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Figure 7. The PDF of the projected flux-decrement p{SJ^a) is plotted. The two different approximations that we have studied produces neai' identical results. Thi'ee 
different redshifts are shown. The lower (higher) redshift corresponds to broader (narrower) distribution. The lines of sight are drawn from a 512^ simulation box. A 
total of 512^ LOS are analysed which are distributed on a regular grid. 

5.2 A Minimal Tree-Model 

The spatial length scales, corresponding to the small angular scales, that are relevant in our discussion, are in the highly non-linear regime. Assuming 
a tree model for the matter correlation hierarchy in the highly non-linear regime, one can write the general form of the A'^th order correlation 
function j'f^ as a product of tw o-point correlation function 5^^' ( Groth & Peebles 1923; iDavis & Peebles 1977 : Frv<1984 iBalian & Schaeffer 19891: 
ISzapudi & Szala^^ll993l,[T997l) In Fourier space such an ansatz means that the hierarchy of multispectra can be written as sums of products of the 
matter power- spectrum: 

Bs{\^^,\^2M)-Zk,^0=Q-i{Ps{k^)P5{k2)+Ps{k2)Ps(k3)+Ps{k'i)P5{k^)) (15) 

Ti(ki,k2,k3,k4)j^fc^=o = Ra P«(A:i)P5(|ki +k2|)P4(|ki +k2 +k3|) +cyc.perm. + Ps{ki)P/,{k2)Ps{k3) + cyc.perm. (16) 

The subscript ki — 0) represents a Dirac's 5 function (^d [X]r=i(^')]- Different hierarchical models differ in the way they specifiy the amplitudes 
Qa, Ra, Rb etc. The results that we will derive are valid for generic hierarchical models. Using Eq.dlSI) and Eq.J16l> we can write: 

(STl), = i3Q3)C3[iXir)f] = S^{STl)l {5Ft)c = (12i?. + ARb)Ci[{X {r)f] = S2{5TI)1 (17) 

Where we have defined the following quantities: 

[iJoir))n - £ ^^^iXir))^^-''dr; ^r)] ^ f ^P. (^) . (18) 

These lower order moments are next used to construct the void probability function (VPF) which also acts as a generating function for the 
normalised Sp parameters for flux-decrement 5iFa{^)- The probability distribution function (PDF) can be constructed from the knowledge of the 
VPF as will be detailed in the next section. 

5.3 Probability Distribution Function of Transmitted Flux 

Before we discuss our techniques to construct the PDF for the projected transmitted flux few explanations regarding the notation and its link to the 
entire hierarchical paradigm are in order. Although our primary aim is to construct the PDF it is usually easier to construct the moment generating 
function ^(y) that will be introduced later. The lower order normalised cumulants (also known as the Sp parameters) up to an arbitrary order can 
be recovered, if necessary, using a Taylor expansion of 4'{y). The function h{x) is related to (jii^x) through an inverse Laplace transform, i.e., a 
knowledge of Sp parameters to an arbitrary order can be used to construct the function h{x). The function h{x) on the other hand is simply a scaled 
PDF and the real PDF can be extracted from it using the definition of scaling variable x (to be defined later). 

The construction of PDF for 5 J- a involves additional steps. First the generating function ^"(y) for the lower order normalised cumulants 5*^ 
is expressed in terms of the underlying (t){y). Next we define a new scaling variable for SJ^a which is denoted by rj. We show that, under certain 
simplifying approximation, the cumulant generating function for 77 i.e. ( j/) is same as the underlying cumulant generating function 4'{y). Once this 
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proved, it is easy see tliat tlie scaled PDFs H^(y) and h(y) too are same as they are simply the Laplace trans form of (y) and (j> ( y) respectively. 
This technique we use has already b een used to derive similar relations in the context of weak lensing statistics jMunshi & Jainl200ol . l200ll : lvalagea3 
bOOdlBemardeau & Valageasll200d) . 

The generating function <I>" (y) (or the Void Probability Function, i.e. VPF) for the lower order normalised cumulants Sp for the flux-decrement 
can be expressed as: 



QCt 



(19) 



Our aim is to relate the VPF ^"{y) and (jjiu) where (f>"{y) is the VPF of the underlying matter distribution. The exact functional form of (jyly) is 
fixed by choosing a specific hierarchical ansatz. We will keep our discussion here completely general without any reference to any specific form of 
the correlation hierarchy. Using Eq.l|17t we can write: 



Using the VPF of the underlying mass distribution 



p=i 



■2\(p-l) 



c T-min 



[x), the above expression can be rewritten in a more compact form: 



dr d 



(20) 



(21) 



We will define a variable rj that will make the analysis simpler. The VPF for 77, denoted "I>''(y), can now be expressed in terms of the VPF of the 
underlying mass distribution (t>{y) by the following expression: 



dr 



\5n 



(22) 



The scaling properties of the PDF in an hierarchical model are encoded in a scaling function. The scaling function (x) is related to the PDF of rj 
i.e.p{ri) where as the scaling function h{x) is associated with the PDF of the underlying mass distribution 5. The following relation relates H^{x) 
and the scaling function h(x). 



)dr 



d\{r){5Fl)c 



h 



x{5Tl).d\{r) 



«;.(r)|5^™|J-(r(r) 



(23) 



The a bove relation is derived using the following definition of (x) in terms of the moment generating function or VPF (y) jSalian & SchaeffeJ 
Il989h as the moment generating function ^''(j/) and the scaled PDF H^{x) are linked through an inverse Laplace transform: 



^exp(xy)$''(y) 



(24) 



A similar result holds for h{x) and <j){y)- Using an approximate form for the integrals we recover the following relation: 

Pc(5T^)=piv)/\SJ'Tl (25) 

In terms of the cumulants this means Sp = S'p/|5J^™'"j''^^. It is instructive to note that the final result is independent of any detailed modelling. 
We will use a lognormal PDF as well as the PDF constructed using the hierarchical ansatz to model p{ri). 

It is important to note that the function h{x) is f undamental to all scsJing analysis. It was initially introduced to study the clustering of galaxies 
later was used in diverse cosmological studies (e.g. IValageas, Schaeffer & Silkl j 19991) ): e.g. the mass function of collapsed objects that includes 
clusters, groups and galaxies can all be described the scaling function h{x). The multiplicity function rj{M, z) at a redshift x for collapsed objects of 
mass M can be expressed via: r]{M, z)dM/M = {p/M)xh{x)dx. Here the scaling function x can be expressed as 2; = [1 + 5{AI, z)]/^2 and p is 
the mean physical density of the Universe. Thus, scaling arguments can also provide an alternat ive to the usual mass function calculation based on the 
Press-Schecter formalism. These issues have been discussed extensively in the literature (e.g. IValageas & SchaeffeJ f20o5h ). The volume averaged 
two-point correlation function ^2 is evaluated at a corresponding relevant scale. It is interesting to notice that the functional form for the scaling 
function h{x) depends only on the form for initial power spectrum and its qualitative behaviour can be predicted using very general arguments. Two 
different asymptotes are particularly well understood; (1) x ^ 1 : h{x) oc x"^^ and (2) a; ^ 1 : h{x) oc x"^^^ exp(— x/x*) with lj ~ 0.5 and 
LJs ~ —3/2. More detailed discussions about some of the salient features are reported in Appendix-A. 

The numerical results for the 2D PDF of rj are presented in FigureQ. The results are depicted for three redshifts, z = 2,3, 4. These results 
were computed using a grid of 512^ LOS drawn from a 512^ simulation box.The underlying cosmology is same as the ones adopted for 3D studies. 
The resulting PDFs computed using this grid are stable down to O(10~*). Notice that the projected PDF is closer to a Gaussian in form than its 3D 
counterpart. 
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Figure 8. The cumulative bias 6(> 5^-"^) of the flux is plotted as a function of the flux. The log-normal distribution and the hierarchical ansatz generates near identical 
results. For a given positive &Fa. the bias is typically higher for a higher redshift. The bias 6(> &Ta) samples the underlying 3D bias associated with 5 i.e. b((5). 



5.4 Joint PDF of Flux-Decrements of neighbouring Quasar Spectra 

The 2PDF is constructed from the cumulant coffelators. The cumulant correlators can be constructed by correlating flux-decrements of neighbouring 
lines of sight and A2. The cumulant correlators are two-point statistics and generalises the usual two-point correlation function to higher order 
that have the ability to probe non-Gaussianity. More interestingly they are related to the bias associated with high flux-decrement regions. 

C^, = ((5J^.(Jlir<5J^a(n2)«),/(5J^^(n))?+«-'(5J^„(ni)5J^4fi2))c. (26) 

Notice that, by construction, Cn = 1 which is due to the particular normalization that is generally adopted for Cpg. At the lower order the cumulant 
correlators can be expressed in terms of the hierarchical amplitudes (i.e. Qs, Ra, Rb, Sa, St, Sc)'- 

{STl{(li)STc,{(l2))c = 2Q3C3lJo"Jr2] = C^jCalXJSi] = CSi{5Tl)c{STc.iCli)STo.iCl2))c, (27) 

{STl{ni)ST{n2)). = {3Ra + 6Rb)C4[{Jo"fjr2] = C^rCi[{Xf J^^] = CS^{5TI)1{T((i^)H^2))., (28) 

{Fl((l,)Fl{(l2))c = ^RbC4{XfJi2] = Cl^2C4{XfJ?2] = CUJ'^)l{5F^{(l{)5F^{(l2))., (29) 

{5Fi{Q.^)STc,{(l2))a = (24Sa + 3656 + ^S^)C5[{Xf J^^] = Cl^CsiXX] = C2i{5Tl)l{ST„{(li)5Tc,{(l2))c, (30) 

{5Tl{(l^)5Tl{(l2))c = (12Sa + QSb)C5[(XfJi2] = CV,^C,[{Xf X] = C^2{STl)l{ST^ini)5T^{{l2))c. (31) 

The expressions for the cumulant correlators given above apply in the large separation limit, i.e. {STa{^i)SJ-a(02)) ^ (STafl)). Though higher- 
order correction terms can be computed using the framework of the hierarchical ansatz, in pr actice the large separation limi t is reached very quickly 
even for nearby lines of sight. A detailed discussion of the tree-amplitudes can be found in dMunshi. Melott. Colesll 19993) . In general the cmulant 
correlators can be expressed as: 

{ST^^{ni)5n{^2))c = c;,C,+,[{Xir)r^''~^X2{r)] = C^,{STlr+'-^{ST^{{l,)ST^{{l2))c. (32) 

Using these expressions we can construct the joint two-point PDF or the bias. To do so we introduce the generating function /3i^'(i/i, ^2) for the 
normalised cumulant correlators: 

«(2)r„ „ ^ ^ yMr'^ _ l yW2 {5T.{Clir3T.,{Cl2r)c 

t^^'^' i-,p^-<iHsnr.-''~Us^c.{a.)5T^{n2)).' 

Using the expression Eq.l|32|l derived using the hierarchical ansatz we can make further progress: 
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Next we introduce the generating function jS^ {y) for the normalised cumulant correlators. While ' (yi , j/2) is associated with Cpq the generating 

yql. 

-2 \ 



function (y) is the generating function for Cpi In 3D the hierarchical ansatz ensures the fcatorization property Cpq — CpiCq 



(35) 



Joir)\5rr\]. (36) 



In terms of the reduced flux i] the generating function for the joint cumulant correlator can be written as: 

These equations are generic and exact. They only depend on the hierarchical ansatz for the higher order correlation functions. To simplify further we 
replace the relevant integrals by the following approximate results: 

- 1 ujUrc)\ f T, ( I \1 ci2 Ia ^Urc)\ f -n t ^ \ v ^ a {\ A rin\ 

~ TrAr / -——Vs[ , , , ) ; £,„ ~ o^r- ,3 , - / 7^nN9P«(3-73) ^^pM • fi'w] ; Ar = r, - r,. (37) 



2 d\{r,)VJ (2^)2 "^d^Crc)']' 2 d\[r,)ij (2^)2 "^d^W' 

This leads to the following expressions which also preserves the factorization properties for the generating functions: 

/3^"(yi,y2) =/3<''(yi)/3<''(y2). (38) 

From the generating function (y) it is possible to construct the scaling function br, (x) that encodes all the in formation about the bias b{ri). The 
functions 6,, (a;) can be recovered using the following integrals in the complex y plane tealian & Schaeffe3l 19891) : 



h^{x)h,-,(x) = dyl3y{y)exp{xy)- b,j{> x)h,,{> x) = dy exp(a;y). (39) 

In addition to the usual scaling functions h{x) and b{x) that are associated with the scaling properties of PDFp(?7) and bias 6(77) we also include their 
cumulative versions h{> x) and b{> x) which are linked to p{> S) and 6(> 5) i.e. the PDF and bias beyond a threshold S. Using the hierarchical 
ansatz allows us to write the 2PDF using the following factorized form: 

p{Vi, V2) = P(??i)p(??2)(l + b{rii)S,l2b{ri2))driidri2. (40) 
Finally, using the definition of rj we can also write the 2PDF of the flux decrement STa for two-different neighbouring lines of sight: 

Pc,iST4ni),sj^c.{fi2)) = p4SJ^c.{^i))p4SJ^c.{n2)){i + b^STc^i^i))^'^^^^ (4i) 

The bias ba for flux-decrement and the bias associated with the underlying mass distribution rj are related by the following expression: 

b^{ST^) = b„{ri)/\STT"\- (42) 

This approximate result also means that we can write C^^ = C^i/|(5J^™'" l''"^ which implies that C^q = C^q/\5J^'^'"\P+''-'^ . Together, Eq.([25) 
and Eg. (142b are sufficient to model the one and two-point distributions of transmitted flux. The fact that the final results are completely independent 
of specific models for the PDF of underlying density contrast indicates their general validity; we can use both the hierarchical model as well as the 
lognormal distribution as a model for the underlying density contrast (5(x) — r]{x.) — 1. 

Though we have only provided the joint PDF of two lines-of-sight directions in this article the formalism can be extended to tackle arbitrary 
number of line of sights. The statistics derived using multiple lines of sight are directly linked to the statistics of hot-spots in flux-decrement maps. 
It is interesting to notice that the bias b{x) introduced here has also been used in previous studies that dealt with clustering of Lyman-a absorbers. 
In a more general context the correlation function of two over-dense objects ^12 can be expressed in terms of the underlying correlation function 
as: (,i2{ri,r2;xi,X2) ~ b{xi)b(x2)^i2iri , r2). Where scaling paramters xi and X2 are the scaling parameters associated with collapsed objects. 
This is a generic outcome of any hierarchical formalism and have been used successfully to compute the bias associated with collapsed objects 
from clusters and groups to galactic size halos ( Valageas & Schaeffer 2000). As was the case for the scaling function h{x) we quote two asymptotes 
for the bias function b(x) given by teemarde au & Schaeffer 1992); (1) for a; < 1; b{x) oc a;'^-"'''^ (2) for a; > 1; b{x) (X X. A more detailed 
discussion is given in Appendix-A. The bias function b{x) is independent of scale, as was the case for the scaled PDF h{x). Information regarding 
a specific scale can be extracted by evaluating the scaling variable x for that scale. 

We note here that lognormal model (see Appendix-lBl) does not belong to the class of other scale-invariant or scaling models. However the 
joint PDF can also be accurately computed using the lognormal model, and in the large separation limit the lognormal model to can be employed to 
compute the integrated bias. 

The computation of the 2PDF or the equivalently the bias was done using a 512^ grid. We computed the number of LOSs where the Flux 
decrement crosses a given threshold which gives us the quantity p{Si,S2)dSidS2. We also use the fractions of LOSs where SJ^a crosses a 

threshold to compute J^°° p{S)d5. Next using the expression given in Eq.Q we compute the resulting bias which is presented in Figure ([8}. 
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6 CONCLUSIONS 

The diffuse intergalactic medium acts as a significant reservoir of baryons at low to intermediate redshift (z < 5), wliicli can be probed via 
the absorption lines in the spectrum of distant QSOs. The Lyman-a absorption lines along the line of sights are due to a vast range of completely 
different classes of object. These objects include under-dense neutral hydrogen clouds, halos of large and over-dense systems which may have already 
reached virial equilibrium, as well as UV heated systems that are strongly coupled to their environments. Previous studies have mainly dealt with 
the problem of detailed modelling of number of these objects as a function of their clustering and internal properties. Using self-consistent scaling 
models, which have a long history and were initially employed in galaxy clustering statistics, several authors have computed the column density 
distribution of Lyman-a absorption systems not only for the low column density Lyman-a forest systems but also for Lyman-limit systems and the 
damped absorption systems. In addition to the hierarchical modelling, lognormal approximation is also quite successful in reproducing the clustering 
statistics of Lyman-a absorption lines. In this study we showed the approach taken by the lognormal approximation and the hierarchical ansatz 
generates near identical results. While previous studies focussed only on one-point statistics we extend these results to the two-point distributions 
and their lower order moments. 

In addition to the statistics of absorption lines, the statistics of the transmitted flux play an important complementary role in Lyman-a studies. 
The flux can be treated as a continuous field. Various statistics which are often employed in analysing the flux include the LOS power spectrum 
measurements, estimation of bispectrum and more recently the entire PDF. The flux PDF contains the information regarding cumulants to an arbitrary 
order. Clearly the flux statistics and the column density distributions are related. One major goal in this study was to unify these two pictures in the 
context of the hierarchical ansatz or scaling models as well as the lognormal approximation. 

The hierarchical model is primarily valid at smaller scales where the correlation functions assumes a hierarchical form. Gravitational clustering 
is known to develop such a form of hierarchy both in the perturbative (quasi-linear) as well as in the highly non-linear regime. The hierarchical ansatz 
can be used to describe the mass functions and the bias associated with collapsed objects using the scaling function h(x) and b{x); the variable a; is a 
scaling variable. Previous studies that have employed the hierarchical ansatz have shown that the column density distribution of absorption systems 
can be described using the same functional form for h{x) and b{x) that are often used for wide ranging studies from the galaxy clustering to thermal 
Sunyaev-Zel'dovich effects or X luminosity of clusters of galaxies. We showed that the same analytical framework can also be used to understand 
the statistics of QSO flux measurements, thereby providing a unified statistical approach to what at first sight appear to be very different observables. 

We have approached the modelling of the flux PDF in two different ways. For the 3D analysis, the PDF of the neutral hydrogen density contrast 
S is modelled according to the hierarchical ansatz as well as lognormal distribution. We find very good agreement with the numerical simulations 
for the entire range of redshifts that we have considered. The predictions from both these models are almost identical and differ only marginally 
in the les s interesting under-dense r egions. We also employ a modified version of scale invariant hierarchical approximation which was developed 
recently ( IValageas & Munshi 2004). This particular approximation provides a very accurate model for the clustering of S. Using the fluctuating 
Gunn-Peterson approximation we map the 5 PDF to that of the transmitted flux Fa- We find reasonable approximation for the allowed range of 
parameters A{z) and /3 that define the Gunn-Peterson approximation. Next we consider the projected or the 2D distribution of flux. For projected 
statistics, we start with linking the cumulants and cumulant correlators for the flux and the underlying neutral hydrogen density contrast S. Then 5 
is taken to be a tracer of the underlying density contrast S, modelled statistically using the hierarchical ansatz. Next it is shown that, under certain 
simplifying assumptions, the PDF of 5 and PDF of a suitably defined reduced flux decrement are linked through a very simple relation. Tests against 
numerical simulations shows good agreement. Results were obtained for the entire bias functions that act as a generating function for the cumulant 
correlators. The scaling function h{x) and the bias function b{x) are also known to describe the number density and bias of over-dense objects 
respectively. 

The formalism developed here can also be used to probe the higher order cross-correlation statistics involving Lyman-a flux decrement and 
weak lensing convergence of CMB maps or those extracted from convergence maps constructed using the weak gravitational lensing of optical 
galaxies. 
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APPENDIX A: HIERARCHICAL ANSATZ: A VERY BRIEF REVIEW 

There is no compl ete analytical model for the evolution of gravitational clustering in the nonlinear regime and the Eulerian and Lagrangian per- 
turbation theories jSemardeau et al|[2002h are often used to describe clustering in the quasi-linear regime. The halo model is extensively used 
and is very successful in modeling fully evolved structure formation. A parallel approach, depends on the hierarchical nature of the correlation 
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functions that develop during tlie nonlinear regime.It has a long history in describing clustering dark matter distribution as well as the collapsed 
objects. Additional assumptions regarding virialization and hydrodynamical equi librium can be used in ass ociation to make specific predictions 
regarding dive rse cosmo logical phenomenon from Luminosity of X-ray cl usters dValaseas & SchaefelbOOOl) . column density distribution of neu- 
tral hydrogen jValageas . Schaeffer & Silk 2002), redshift evolution bias dValageas. Silk & SchaeffeJ ] 1200 ih . mass and luminosity distribution of 
galaxies and clusters fValageas & Schaeffer 1999) , cosmic m icrowave background (CMB) secondaries such as the thermal and kinetic Sunyaev 
and Zel'dovic h effect l,Munshi et a l. 2011; Valag eas. Balbi & S ilk 2001), reheating and reionization of th e Universe (|Valageas & Schaeffer 199^; 
Valageas. Sc haeffer & Silk|l2002l) . Weak lensing observables h ave already been studied in this framework ( IValageasll200ll ; Bemardeau & Valageaj 
200Q : iMunshi & Jainll200ot 1200 iV Valageas & Schaeffeilll999l). 

The hierarchical ansatzTBalian & Schaeffer 1989) is remarkably successful at making app roximate calculations o f the entire PDF and bias of 
the density field, improving significantly upon the order-by-order analysis of other approaches jBemardeau et all ( |2002|) for a detailed review). The 
hierarchical ansatz in the highly non-linear regime depends on assuming a specific correlation hierarchy, on the other hand in the quasi-linear regime 
it can be linked to the gravitational dynamics using perturbative analysis. 



Al Highly Non-linear Regime 

The PDFp(5) and the bias b{5) can both be constructed from the knowledge of the VPF (/>(t/) = Spy^ /p\ and its two-point analog r(t/) = 

X]p Cpiy^ /p\. Where the parameters Sp and Cpi are normalized cumulants and cumulant correlators for the density field. The PDF is related to the 
scaling function h(x) in the text of the paper p{5) = h{x)/^ with x — {1 + 5) /^2- 

dy I {1 + 5)y - (t>{y) 



i^exp[i 



6 



]; 6(5)p(5)= ^r(y)exp[ 



- (1 + 5)y - <l>{y) ^ 



(Al) 



The modelling of (j}{y) and T{y) needs a detailed knowledge of the entire correlation hierarchy. The detailed knowledge of the entire correlation 
hierarchy is encoded in the vertex generating function G{t)). Typically for large values of y the VPF exhibits a power law (j)(y) = ay^^'^. There 
are no theoretical estimates of cj and it is generally estimated from numerical simulations. The parameter typically takes a value uj = 0.3 for CDM 
like spectra. For small but negative values of y the functions (j>{y) ™d '''(v) develops a singularity in the complex plane which is described by the 
following parametrization. 

4>{y) ^ (l^s - asr{ujs){y - ysy^'; T{y) ^ Ts ~ bs{y - ysy^"~^ ■ (A2) 
The singularity structure of (f>iy) ™d ''"(j/) depends on the nature of the vertex generating function G(r) and its behaviour near the singularity Ts'. 



r(-i/2) 



g'irs)g"irs) 



g"'{rs) 



3/2 



2g'irs)g"{r, 



g"'irs) 



1/2 



(A3) 



On the other hand the parameters cj and a can be described in terms of a parameter ys which in turn describes the exponential decay of the PDF for 
large density contrast S: 



OJ = ka/{ka+2), 



-k: 



fea/fca+2. 



'ys A:, + 



(A4) 



The PDF and the bias thus has two distinct regimes that are dictated by the two asymptotes. For intermediate values of 5 the PDF shows a power law 
behaviour. The PDF and the bias are given by the following expression: 



<<l + 5<<,t; p{5) 



l-uj fl + 5 



1/2 



r(a 



1 + 5 



For large values of 5 the PDF on the other hand shows an exponential behaviour: 



C2 ^ 6 ^ ^ 2;*t2 ^ 



(l-")/2 



g'{rs) 6 

At very small values of 5 the PDF shows an exponential decay which is described only by the parameter uj: 



1 1 x ^ ^ c ix\ -i/(i-ij)ftj/{i-uj) / (1— o;)!/" r 



1-UJ 



1/2 



1 -CJ 



(l-c»j)/2aj 



(A5) 



(A6) 



(A7) 



The range of S for which the power law regime is valid depends on the value of ^2. For smaller values of ^2 the power law regime is less pronounced. 
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A2 Quasi-linear Regime 

In the quasi-linear regime the parameters that describe the generating functions can be directly Iniked with gravitational dynamics teemardeaulll99l 
I1994I) . The PDF and bias in the intermediate power-law regime in the quasi-linear regime is given by the following result: 

Qi>{T)^g{T)^l = 5. (A9) 
For the large density contrast the exponential decay take the following form: 

p(5)d,.!^,-/^exp[-|y.|l + ^ld5; 6(5) = - 1(1 + ^. (AlO) 
4V7i" L 42 ^2 J y'{rs) 42 

These results ignore the loop corrections to tree-level perturbation theory but take into account correlation hierarchy to an arbitrary order. 



APPENDIX B: LOGNORMAL DISTRIBUTION 



The evolution of the P DF p(5) of density field 5 due to gravitational clustering has been studied ex tensively in many cosmological context. The 
lognormal distribution faamiltoiJ[T98^ : IColes & Joneslll99ll : | Bouchet et ai l 1 9931 ; IXofman et al.|[T994l) is known to accurately reproduce the results 
from numerical simulation in the quasi-linear regime and has been used to model both one- and two-point PDF. In general the lognormal distribution 
can provide a good statis tical description of a random variable if it can be modeled as a product of many other random variables tKavo. Taruva. Sutol 
l200ll : lTaruva et alj2002l) . 



p{5)dS = —r= exp 



a! 

'2E 



dS 
1 + 5' 



/27rE 

E = ln(l + cr^); A = ln[(l + 5)7(1+72); 



p{Si, S2)dSidS2 
A, = ln[(l + 5,)V(l + 'T2); 



^12 



exp 



E(A? + Ai) - 2A:i2AiA2 



2(E2 



X?2) 



dSi dS2 
1 + 5i 1 + 52 ' 



A:i2 = ln(l + ^12) 

In the lowest order in ^12 the lognormal model takes a factorized form. In this limit a bias function 6(5) can be defined. 

p(5i, 52) = p(5i)p(52)[l + b{5i)^i2b{52)]; 6(5,) = A,/E. 
The lognormal distribution has also been used to model many cosmological observation from galaxy surveys to weak lensing data. 
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